Este trabajo propone una metodología innovadora para construir redes de dependencia entre acciones en el mercado financiero mediante el uso de Graphical LASSO (GLASSO). Al modelar las acciones como nodos de una red y analizar sus precios a lo largo del tiempo, buscamos identificar patrones de dependencias condicionales que podrían revelar tendencias y anomalías de interés para la toma de decisiones de inversión. La propuesta se centra en proporcionar una herramienta visual para comprender la estructura de dependencia en los activos financieros.
required.packages <- c('glasso', 'colorRamps', 'igraph', 'RColorBrewer', 'threejs', 'htmlwidgets','quantmod', 'tidyverse', 'PerformanceAnalytics')
new.packages <- required.packages[!(required.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) install.packages(new.packages, repos='http://cran.us.r-project.org')
library(ggplot2)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(glasso)
## Warning: package 'glasso' was built under R version 4.3.1
library(colorRamps)
## Warning: package 'colorRamps' was built under R version 4.3.3
library(igraph)
## Warning: package 'igraph' was built under R version 4.3.3
##
## Attaching package: 'igraph'
## The following objects are masked from 'package:dplyr':
##
## as_data_frame, groups, union
## The following objects are masked from 'package:stats':
##
## decompose, spectrum
## The following object is masked from 'package:base':
##
## union
library(RColorBrewer)
library(threejs)
## Warning: package 'threejs' was built under R version 4.3.3
library(htmlwidgets)
library(quantmod)
## Warning: package 'quantmod' was built under R version 4.3.3
## Loading required package: xts
## Warning: package 'xts' was built under R version 4.3.3
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
## ######################### Warning from 'xts' package ##########################
## # #
## # The dplyr lag() function breaks how base R's lag() function is supposed to #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or #
## # source() into this session won't work correctly. #
## # #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop #
## # dplyr from breaking base R's lag() function. #
## # #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning. #
## # #
## ###############################################################################
##
## Attaching package: 'xts'
## The following objects are masked from 'package:dplyr':
##
## first, last
## Loading required package: TTR
## Warning: package 'TTR' was built under R version 4.3.3
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ lubridate 1.9.2 ✔ tibble 3.2.1
## ✔ purrr 1.0.1 ✔ tidyr 1.3.0
## ✔ readr 2.1.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ lubridate::%--%() masks igraph::%--%()
## ✖ tibble::as_data_frame() masks igraph::as_data_frame(), dplyr::as_data_frame()
## ✖ purrr::compose() masks igraph::compose()
## ✖ tidyr::crossing() masks igraph::crossing()
## ✖ dplyr::filter() masks stats::filter()
## ✖ xts::first() masks dplyr::first()
## ✖ dplyr::lag() masks stats::lag()
## ✖ xts::last() masks dplyr::last()
## ✖ purrr::simplify() masks igraph::simplify()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(PerformanceAnalytics)
## Warning: package 'PerformanceAnalytics' was built under R version 4.3.3
##
## Attaching package: 'PerformanceAnalytics'
##
## The following object is masked from 'package:graphics':
##
## legend
Se analizaran cerca de ~200 stocks, tomando algunos de los más conocidos a nivel internacional:
symbols <- c("ABT", "ACN", "ATVI", "ADBE", "AFL", "A", "APD", "ALL", "GOOG", "MO", "AMZN",
"AEP", "AXP", "AIG", "AMT", "AMP", "AMGN", "ADI", "AON", "AAPL", "APTV", "T",
"ADSK", "ADP", "AZO", "BAX", "BDX", "BIIB", "BLK", "BMY", "AVGO", "CDNS",
"CAT", "CNC", "SCHW", "CVX", "CB", "CI", "CTAS", "CSCO", "C", "CME", "KO",
"CTSH", "CL", "CMCSA", "ED", "STZ", "GLW", "COST", "CCI", "CVS", "DHI",
"DHR", "DE", "DLR", "DG", "DLTR", "D", "DTE", "DUK", "ETN", "EBAY", "ECL",
"EW", "EMR", "EOG", "EQR", "ES", "EXC", "XOM", "FAST", "FIS", "FE", "FISV",
"GE", "GM", "GILD", "GS", "HCA", "HLT", "HD", "HON", "HPQ", "HUM", "IDXX",
"ITW", "ILMN", "INTC", "ICE", "IBM", "INTU", "ISRG", "JNJ", "JPM", "KMB",
"KLAC", "KHC", "KR", "LRCX", "LLY", "LMT", "LOW", "MPC", "MAR", "MMC", "MA",
"MCD", "MDT", "MRK", "MET", "MTD", "MCHP", "MSFT", "MDLZ", "MNST", "MCO",
"MS", "MSI", "NFLX", "NEM", "NEE", "NKE", "NSC", "NOC", "NVDA", "ORLY", "OXY",
"ORCL", "PCAR", "PAYX", "PYPL", "PEP", "PFE", "PM", "PSX", "PNC", "PPG", "PG",
"PGR", "PLD", "PRU", "PEG", "PSA", "QCOM", "REGN", "ROP", "ROST", "SPGI",
"SBAC", "SRE", "SHW", "SPG", "SO", "SBUX", "SYK", "SYY", "TGT", "TEL", "TXN",
"BK", "TRV", "DIS", "TMO", "TJX", "TSN", "USB", "UNP", "UNH", "UPS", "VLO",
"VRSK", "VZ", "V", "WMT", "WBA", "WM", "WEC", "WFC", "WMB", "YUM", "ZTS")
sectors <- c("Health Care", "Information Technology", "Information Technology",
"Information Technology", "Financials", "Health Care", "Materials",
"Financials", "Information Technology", "Consumer Staples", "Consumer Discretionary",
"Utilities", "Financials", "Financials", "Real Estate", "Financials",
"Health Care", "Information Technology", "Financials", "Information Technology",
"Consumer Discretionary", "Telecommunication Services", "Information Technology",
"Information Technology", "Consumer Discretionary", "Health Care", "Health Care",
"Health Care", "Financials", "Health Care", "Information Technology",
"Information Technology", "Industrials", "Health Care", "Financials",
"Energy", "Financials", "Health Care", "Industrials", "Information Technology",
"Financials", "Financials", "Consumer Staples", "Information Technology",
"Consumer Staples", "Consumer Discretionary", "Utilities", "Consumer Staples",
"Information Technology", "Consumer Staples", "Real Estate", "Consumer Staples",
"Consumer Discretionary", "Health Care", "Industrials", "Real Estate",
"Consumer Discretionary", "Consumer Discretionary", "Utilities", "Utilities",
"Utilities", "Industrials", "Information Technology", "Materials",
"Health Care", "Industrials", "Energy", "Real Estate", "Utilities",
"Utilities", "Energy", "Industrials", "Information Technology",
"Utilities", "Information Technology", "Industrials", "Consumer Discretionary",
"Health Care", "Financials", "Health Care", "Consumer Discretionary",
"Consumer Discretionary", "Industrials", "Information Technology",
"Health Care", "Health Care", "Industrials", "Health Care", "Information Technology",
"Financials", "Information Technology", "Information Technology",
"Health Care", "Health Care", "Financials", "Consumer Staples",
"Information Technology", "Consumer Staples", "Consumer Staples",
"Information Technology", "Health Care", "Industrials", "Consumer Discretionary",
"Energy", "Consumer Discretionary", "Financials", "Information Technology",
"Consumer Discretionary", "Health Care", "Health Care", "Financials",
"Health Care", "Information Technology", "Information Technology",
"Consumer Staples", "Consumer Staples", "Financials", "Financials",
"Information Technology", "Information Technology", "Materials",
"Utilities", "Consumer Discretionary", "Industrials", "Industrials",
"Information Technology", "Consumer Discretionary", "Energy",
"Information Technology", "Industrials", "Information Technology",
"Information Technology", "Consumer Staples", "Health Care",
"Consumer Staples", "Energy", "Financials", "Materials", "Consumer Staples",
"Financials", "Real Estate", "Financials", "Utilities", "Real Estate",
"Information Technology", "Health Care", "Industrials", "Consumer Discretionary",
"Financials", "Real Estate", "Utilities", "Materials", "Real Estate",
"Utilities", "Consumer Discretionary", "Health Care", "Consumer Staples",
"Consumer Discretionary", "Information Technology", "Information Technology",
"Financials", "Financials", "Consumer Discretionary", "Health Care",
"Consumer Discretionary", "Consumer Staples", "Financials", "Industrials",
"Health Care", "Industrials", "Energy", "Industrials", "Telecommunication Services",
"Information Technology", "Consumer Staples", "Consumer Staples",
"Industrials", "Utilities", "Financials", "Energy", "Consumer Discretionary",
"Health Care")
data <- data.frame(symbols = symbols, sectors = sectors)
data
La elección de la ventana de tiempo adecuada para el análisis de clustering de valores bursátiles y sus retornos depende de los objetivos específicos del análisis y de las características de los datos que deseas capturar. Aquí algunos enfoques para elegir una ventana de tiempo prudente:
Para evaluar estos valores en diferentes condiciones de mercado (por ejemplo, tanto en tiempos de expansión como de recesión), usaremos una ventana de 36 meses.
# end_date <- Sys.Date() - days(1) # ayer
# start_date <- end_date - months(36)
# print(end_date)
# print(start_date)
Se define una función que utiliza la función getSymbols() pertenece al paquete quantmod y se utiliza para descargar datos financieros de diferentes fuentes (por defecto Yahoo Finance). Luego, la función Cl(), del mismo paquete, extrae los precios de cierre ajustado de un objeto tipo OHLCV (Open, High, Low, Close, Volume).
get_stock_data <- function(symbol) {
tryCatch({
stock_data <- getSymbols(symbol, src = "yahoo", from = start_date, to = end_date, auto.assign = FALSE)
stock_data <- Cl(stock_data)
colnames(stock_data) <- symbol
return(stock_data)
}, error = function(e) {
message(paste("No se pudieron obtener datos para", symbol))
return(NULL)
})
}
Se aplica la función creada para las stocks definidas previamente, aplicando un merge para combinar las series temporales de cada símbolo en un único xts (serie temporal extensible).
# stock_list <- lapply(symbols, get_stock_data)
# stock_data <- do.call(merge, stock_list)
Para evitar volver a consultar la API de yahoo finance:
# stock_data <- as.data.frame(stock_data)
# stock_data$fecha <- rownames(stock_data)
# write.table(stock_data, file = sprintf("yahoo_data_%s.txt", end_date), row.names = F, col.names = T)
Lectura de los datos guardados con fecha 2024-12-04, con 36 meses de valores:
# stock_data <- read.table(sprintf("yahoo_data_%s.txt", end_date), header = TRUE)
stock_data <- read.table("yahoo_data_2024-12-04.txt", header = TRUE)
# Convertir la columna 'fecha' a nombres de fila
rownames(stock_data) <- stock_data$fecha
stock_data$fecha <- NULL
# Filtrar por stocks bajo estudio
stock_data <- stock_data %>% select(intersect(names(stock_data), symbols))
# Mostrar el resultado
stock_data
Se descartan las columnas con datos faltantes o rotas. Luego, se preparan los dataframes para analizar los datos obtenidos de Internet.
Además, se calculan los retornos usando la función ROC que proviene del paquete TTR, el cual usualmente se instala de manera automática junto con quantmod, y sirve para calcular la tasa de cambio de una serie temporal de precios u otros datos financieros. Su uso típico es con precios de un activo en diferentes momentos y se emplea frecuentemente en análisis técnico para medir cuánto varía el precio entre un punto actual y un punto anterior.
# Eliminar columnas con datos faltantes
# stock_data <- na.omit(stock_data) # quita muchas filas, por ausencias o falta de nombres
# Calcular retornos diarios
returns <- ROC(stock_data, type = "discrete")
# Convertir a tibble y agregar columnas útiles para el análisis de valores
stock_data_df <- as.data.frame(stock_data) %>%
rownames_to_column(var = "Date") %>%
mutate(Date = ymd(Date))
# Convertir a tibble y agregar columnas útiles para el análisis de retornos
returns_df <- as.data.frame(returns) %>%
rownames_to_column(var = "Date") %>%
mutate(Date = ymd(Date))
Un rápido vistazo a los datos:
print(summary(stock_data_df[-1][0:9]))
## AAPL MSFT GOOG AMZN
## Min. :125.0 Min. :214.2 Min. : 83.49 Min. : 81.82
## 1st Qu.:156.0 1st Qu.:271.9 1st Qu.:110.75 1st Qu.:115.01
## Median :173.3 Median :324.9 Median :132.31 Median :140.39
## Mean :176.3 Mean :331.8 Mean :132.32 Mean :142.74
## 3rd Qu.:189.9 3rd Qu.:406.3 3rd Qu.:148.68 3rd Qu.:173.16
## Max. :242.6 Max. :467.6 Max. :192.66 Max. :214.10
## V JNJ WMT JPM
## Min. :177.7 Min. :144.4 Min. :39.43 Min. :102.0
## 1st Qu.:213.4 1st Qu.:156.6 1st Qu.:47.22 1st Qu.:132.2
## Median :231.3 Median :162.8 Median :51.39 Median :145.8
## Mean :238.2 Mean :163.3 Mean :54.65 Mean :157.1
## 3rd Qu.:267.9 3rd Qu.:170.7 3rd Qu.:59.35 3rd Qu.:184.4
## Max. :316.6 Max. :186.0 Max. :93.51 Max. :250.3
## MA
## Min. :283.4
## 1st Qu.:352.3
## Median :381.5
## Mean :395.8
## 3rd Qu.:444.8
## Max. :532.9
stock_data_df
print(summary(returns_df[-1][0:9]))
## AAPL MSFT GOOG
## Min. :-0.0586795 Min. :-0.0771563 Min. :-0.0963500
## 1st Qu.:-0.0083559 1st Qu.:-0.0089237 1st Qu.:-0.0116156
## Median : 0.0008576 Median : 0.0006171 Median : 0.0014792
## Mean : 0.0006590 Mean : 0.0005239 Mean : 0.0004541
## 3rd Qu.: 0.0100771 3rd Qu.: 0.0107313 3rd Qu.: 0.0113377
## Max. : 0.0889746 Max. : 0.0822681 Max. : 0.0996518
## NA's :1 NA's :1 NA's :1
## AMZN V JNJ
## Min. :-0.1404944 Min. :-0.0549146 Min. :-0.0398326
## 1st Qu.:-0.0127029 1st Qu.:-0.0070720 1st Qu.:-0.0056965
## Median : 0.0003105 Median : 0.0013160 Median : 0.0000632
## Mean : 0.0005820 Mean : 0.0006766 Mean :-0.0000361
## 3rd Qu.: 0.0139554 3rd Qu.: 0.0074102 3rd Qu.: 0.0056167
## Max. : 0.1353590 Max. : 0.1059908 Max. : 0.0607282
## NA's :1 NA's :1 NA's :1
## WMT JPM MA
## Min. :-0.113757 Min. :-0.0646779 Min. :-0.0562243
## 1st Qu.:-0.005197 1st Qu.:-0.0075801 1st Qu.:-0.0068125
## Median : 0.001094 Median : 0.0012035 Median : 0.0008128
## Mean : 0.001019 Mean : 0.0006883 Mean : 0.0007148
## 3rd Qu.: 0.007802 3rd Qu.: 0.0089027 3rd Qu.: 0.0083668
## Max. : 0.069865 Max. : 0.1154454 Max. : 0.0912333
## NA's :1 NA's :1 NA's :1
returns_df
Para mostrar el valor de la técnica Graphical LASSO, en este caso de uso, vamos a aplicar un análisis de correlación de los valores.
Tomando solo 30 stocks stocks al azar:
# Seleccionar 30 variables al azar de stock_data_df excluyendo la primera columna
set.seed(42) # Fijar semilla para reproducibilidad
random_columns <- sample(names(stock_data_df)[-1], 33)
# Crear un subconjunto del dataframe con las columnas seleccionadas al azar
stock_data_random <- stock_data_df[random_columns]
# Calcular la matriz de correlación con las 30 columnas seleccionadas
cor_matrix_values <- cor(stock_data_random, use = "pairwise.complete.obs")
# Visualizar la matriz de correlación
corrplot::corrplot(cor_matrix_values, method = "color", tl.cex = 0.6, tl.col = "black")
OBSERVACIÓN: hay muchos valores altamente correlacionados, tanto
negativa como positivamente.
Preparando un filtro para las variables muy correlacionadas (abs > 0.9), para visualizar mejor los datos:
# Filtrar símbolos con alta correlación en valores
pos_corr_pairs_values <- which(cor_matrix_values > 0.8 & abs(cor_matrix_values) < 1, arr.ind = TRUE)
neg_corr_pairs_values <- which(cor_matrix_values < -0.8 & abs(cor_matrix_values) < 1, arr.ind = TRUE)
# Asumimos que ya existen las variables pos_corr_pairs_values y neg_corr_pairs_values
# provenientes del código que mostraste.
# Crear un data frame con los pares altamente correlacionados positivamente
pos_pairs <- data.frame(
Symbol1 = rownames(cor_matrix_values)[pos_corr_pairs_values[, 1]],
Symbol2 = colnames(cor_matrix_values)[pos_corr_pairs_values[, 2]],
Correlation = cor_matrix_values[pos_corr_pairs_values]
)
# Crear un data frame con los pares altamente correlacionados negativamente
neg_pairs <- data.frame(
Symbol1 = rownames(cor_matrix_values)[neg_corr_pairs_values[, 1]],
Symbol2 = colnames(cor_matrix_values)[neg_corr_pairs_values[, 2]],
Correlation = cor_matrix_values[neg_corr_pairs_values]
)
# Si deseas, puedes ver una muestra de estos data frames:
head(pos_pairs)
head(neg_pairs)
Ahora bien, que sectores corresponden a estas correlaciones:
subset(data, symbols %in% c("MMC", "HLT"))
subset(data, symbols %in% c("MMC", "UPS"))
Encontramos casos correlaciones que parecen no tener una colinealidad o relación directa, según su sector.
Filtrando dataframes para un EDA:
# Columnas con correlacion altamente positiva
pos_corr_stocks <- stock_data_df %>% select(Date, MMC, HLT)
# Columnas con correlacion altamente positiva
neg_corr_stocks <- stock_data_df %>% select(Date, MMC, UPS)
stock_data_df
Correlación positiva:
price_trend_plot_filtered <- pos_corr_stocks %>%
pivot_longer(-Date, names_to = "Symbol", values_to = "Close") %>%
ggplot(aes(x = Date, y = Close, color = Symbol)) +
geom_line() +
labs(title = "Algunas tendencia de precios de cierre", x = "Fecha", y = "Precio de Cierre") +
theme_minimal()
print(price_trend_plot_filtered)
Correlación negativa:
price_trend_plot_filtered <- neg_corr_stocks %>%
pivot_longer(-Date, names_to = "Symbol", values_to = "Close") %>%
ggplot(aes(x = Date, y = Close, color = Symbol)) +
geom_line() +
labs(title = "Algunas tendencia de precios de cierre", x = "Fecha", y = "Precio de Cierre") +
theme_minimal()
print(price_trend_plot_filtered)
Correlación positiva:
return_hist_plot_filtered <- pos_corr_stocks %>%
pivot_longer(-Date, names_to = "Symbol", values_to = "Return") %>%
ggplot(aes(x = Return, fill = Symbol)) +
geom_histogram(bins = 50, color = "black", alpha = 0.7) +
labs(title = "Histograma de sus retornos diarios", x = "Retorno diario", y = "Frecuencia") +
theme_minimal()
print(return_hist_plot_filtered)
Correlación negativa:
return_hist_plot_filtered <- neg_corr_stocks %>%
pivot_longer(-Date, names_to = "Symbol", values_to = "Return") %>%
ggplot(aes(x = Return, fill = Symbol)) +
geom_histogram(bins = 50, color = "black", alpha = 0.7) +
labs(title = "Histograma de sus retornos diarios", x = "Retorno diario", y = "Frecuencia") +
theme_minimal()
print(return_hist_plot_filtered)
CONCLUSION: Analizar valores bursátiles/stocks e intentar relacionarlos por sus correlaciones o sus distribuciones de ROC es confuso, sobretodo para el ojo inexperto. Conocer el tipo de industria o negocio que sustenta al stock para relevante, pero aún así no parece suficiente para agruparlos por su compartamiento bursatil.
Aplicación de la técnica, paso por paso, partiendo del dataframe de retornos en lugar de los valores de cierre. Esto, resulta ventajoso para:
cov_matrix <- cov(na.omit(returns)) # al omitir na se pierden varios años de datos VER
La matriz de covarianza contiene la información estadística de cómo las variables (valores o retornos) varian conjuntamente. La matriz de correlación usa a la covarianza y trata de indicar la intensidad y la dirección de la relación entre variables, tal como analizamos anteriormente.
Sin embargo, existen correlaciones espureas o indirectas, indescifrables a partir de una matriz de correlación.
Graphical LASSO aborda este problema trabajando sobre la matriz de precisión que, para cada elemento \(\Omega_{jk}\) describe la dependencia condicional entre dos variables (\(X_j\) y \(X_k\)), es decir, cuánto están relacionadas cuando se eliminan los efectos de todas las demás variables.
Luego, la matriz de precisión \(\Omega\) puede ser aproximada por la inversa de la matriz de covarianza \(\Sigma^-1\). Entonces, Graphical LASSO estima una matriz de precisión esparsa \(\hat{\Omega}\), que es la inversa de la matriz de covarianza regularizada con la siguiente función de optimización:
\[ \hat{\Omega} = \underset{\Omega \succ 0}{\arg \min} \ \Big[ \text{tr}(\Sigma \Omega) - \log\det(\Omega) + \lambda \|\Omega\|_1 \Big] \]
Donde: - \(\Sigma\): matriz de covarianza. - \(\Omega\): matriz de precisión (inversa de la covarianza). - \(\text{tr}(\Sigma \Omega)\): traza del producto entre \(\Sigma\) y \(\Omega\). - \(\log\det(\Omega)\): logaritmo del determinante de \(\Omega\), que garantiza que \(\Omega\) sea positiva definida. - \(\|\Omega\|_1\): suma de los valores absolutos de los elementos de \(\Omega\), que induce esparsidad. - \(\lambda\): parámetro de regularización (también llamado \(\rho\) en algunos contextos).
Intuición del término de regularización
A medida que \(\lambda\) aumenta, más elementos de \(\hat{\Omega}\) se fuerzan a cero, eliminando conexiones entre variables y generando un grafo más esparso. Un \(\lambda\) bajo permite conexiones más débiles y genera un grafo más denso.
Usamos la función glasso para resolver este problema. El
resultado incluye la matriz de precisión regularizada (\(\hat{\Omega}\)):
rho <- 0.1 # Valor de regularización bajo
glasso_result <- glasso(cov_matrix, rho = rho)
# Matriz de precisión
precision_matrix <- glasso_result$wi
La correlación parcial mide la relación directa entre dos variables \(X_j\) y \(X_k\), eliminando los efectos de todas las demás variables (\(X_{\setminus\{j,k\}}\)). Se calcula a partir de la matriz de precisión regularizada (\(\Omega\)).
Dado que \(\Omega\) es la matriz de precisión, la correlación parcial entre \(X_j\) y \(X_k\) está dada por:
\[ R_{jk} = -\frac{\Omega_{jk}}{\sqrt{\Omega_{jj} \cdot \Omega_{kk}}} \]
Donde: - \(\Omega_{jk}\): elemento de la matriz de precisión que representa la dependencia condicional entre \(X_j\) y \(X_k\). - \(\Omega_{jj}\) y \(\Omega_{kk}\): elementos diagonales de \(\Omega\), que representan la varianza condicional de \(X_j\) y \(X_k\).
Intuición matemática
Si \(\Omega_{jk}\) es cero, \(R_{jk} = 0\). Esto indica que \(X_j\) y \(X_k\) son condicionalmente independientes, dado el resto de las variables. Si \(\Omega_{jk} \neq 0\), el signo y la magnitud de \(R_{jk}\) determinan la dirección y fuerza de la relación directa entre \(X_j\) y \(X_k\).
# Calcular matriz de correlación parcial
partial_corr_matrix <- matrix(0, nrow = nrow(precision_matrix), ncol = ncol(precision_matrix))
for (j in 1:nrow(precision_matrix)) {
for (k in 1:ncol(precision_matrix)) {
partial_corr_matrix[j, k] <- -precision_matrix[j, k] / sqrt(precision_matrix[j, j] * precision_matrix[k, k])
}
}
diag(partial_corr_matrix) <- 0 # Evitar autoconexiones
colnames(partial_corr_matrix) <- colnames(returns)
rownames(partial_corr_matrix) <- colnames(returns)
Usamos igraph para construir y visualizar un gráfico basado en las conexiones significativas en la matriz de correlación parcial.
# Construir grafo desde la matriz de correlación parcial
stock_graph <- graph_from_adjacency_matrix(partial_corr_matrix > 0.1, mode = "max", weighted = TRUE)
# Asignar sectores a los nodos
V(stock_graph)$sector <- data$sectors[match(V(stock_graph)$name, data$symbols)]
# Crear un vector de colores únicos para cada sector
unique_sectors <- unique(V(stock_graph)$sector)
sector_colors <- rainbow(length(unique_sectors)) # Genera una paleta de colores
# Mapear cada sector a un color
sector_color_map <- setNames(sector_colors, unique_sectors)
# Asignar colores a los vértices según su sector
V(stock_graph)$color <- sector_color_map[V(stock_graph)$sector]
plot(stock_graph,
vertex.size = 10,
vertex.label = NULL,
edge.width = E(stock_graph)$weight * 2,
vertex.color = V(stock_graph)$color,
main = "Grafo de acciones - según retornos y coloreado sector")
Gráfico HTML interactivo con graphjs:
# make interactive graph
stock_graph_js <- graphjs(g=stock_graph,
layout_with_fr(stock_graph, weights=30*E(stock_graph)$width, dim=3),
vertex.size=0.7,
vertex.frame.color="white",
vertex.frame.width=0.2,
vertex.label=paste(names(V(stock_graph)), "-", V(stock_graph)$sector), # Asignar nombres que incluyan el sector
vertex.color = V(stock_graph)$color,
brush=TRUE, # resalte del nodo con un click
showLabels=TRUE, # mostrar nombres con hover
edge.alpha=0.6, # transparencia de las conexiones
bg="black", # background
main="Stocks graph using returns with Graphical Lasso")
# save graph
graph_filename <- paste0("returns_graph_rho_", rho, "_with_sectors.html")
saveWidget(stock_graph_js, file=graph_filename)
# open in browser
browseURL(graph_filename)
OBSERVACIÓN; con returns el método falla y no encuentra relaciones.
Aplicación de la técnica, paso por paso, partiendo del dataframe de valores de cierre de los stocks.
cov_matrix <- cov(na.omit(stock_data)) # al omitir na se pierden varios años de datos VER
Usamos la función glasso para resolver este problema. El
resultado incluye la matriz de precisión regularizada (\(\hat{\Omega}\)):
rho <- 0.9 # Valor de regularización elevado
glasso_result <- glasso(cov_matrix, rho = rho)
# Matriz de precisión
precision_matrix <- glasso_result$wi
La correlación parcial mide la relación directa entre dos variables \(X_j\) y \(X_k\), eliminando los efectos de todas las demás variables (\(X_{\setminus\{j,k\}}\)). Se calcula a partir de la matriz de precisión regularizada (\(\Omega\)).
Dado que \(\Omega\) es la matriz de precisión, la correlación parcial entre \(X_j\) y \(X_k\) está dada por:
\[ R_{jk} = -\frac{\Omega_{jk}}{\sqrt{\Omega_{jj} \cdot \Omega_{kk}}} \]
Donde: - \(\Omega_{jk}\): elemento de la matriz de precisión que representa la dependencia condicional entre \(X_j\) y \(X_k\). - \(\Omega_{jj}\) y \(\Omega_{kk}\): elementos diagonales de \(\Omega\), que representan la varianza condicional de \(X_j\) y \(X_k\).
Intuición matemática
Si \(\Omega_{jk}\) es cero, \(R_{jk} = 0\). Esto indica que \(X_j\) y \(X_k\) son condicionalmente independientes, dado el resto de las variables. Si \(\Omega_{jk} \neq 0\), el signo y la magnitud de \(R_{jk}\) determinan la dirección y fuerza de la relación directa entre \(X_j\) y \(X_k\).
# Calcular matriz de correlación parcial
partial_corr_matrix <- matrix(0, nrow = nrow(precision_matrix), ncol = ncol(precision_matrix))
for (j in 1:nrow(precision_matrix)) {
for (k in 1:ncol(precision_matrix)) {
partial_corr_matrix[j, k] <- -precision_matrix[j, k] / sqrt(precision_matrix[j, j] * precision_matrix[k, k])
}
}
diag(partial_corr_matrix) <- 0 # Evitar autoconexiones
colnames(partial_corr_matrix) <- colnames(returns)
rownames(partial_corr_matrix) <- colnames(returns)
Usamos igraph para construir y visualizar un gráfico basado en las conexiones significativas en la matriz de correlación parcial:
# Construir grafo desde la matriz de correlación parcial
stock_graph <- graph_from_adjacency_matrix(partial_corr_matrix > 0.1, mode = "max", weighted = TRUE)
# Asignar sectores a los nodos
V(stock_graph)$sector <- data$sectors[match(V(stock_graph)$name, data$symbols)]
# Crear un vector de colores únicos para cada sector
unique_sectors <- unique(V(stock_graph)$sector)
sector_colors <- rainbow(length(unique_sectors)) # Genera una paleta de colores
# Mapear cada sector a un color
sector_color_map <- setNames(sector_colors, unique_sectors)
# Asignar colores a los vértices según su sector
V(stock_graph)$color <- sector_color_map[V(stock_graph)$sector]
plot(stock_graph,
vertex.size = 10,
vertex.label = NULL,
edge.width = E(stock_graph)$weight * 2,
vertex.color = V(stock_graph)$color,
main = "Grafo de acciones - según valores y coloreado sector")
Gráfico HTML interactivo con graphjs:
# make interactive graph
stock_graph_js <- graphjs(g=stock_graph,
layout_with_fr(stock_graph, weights=30*E(stock_graph)$width, dim=3),
vertex.size=0.7,
vertex.frame.color="white",
vertex.frame.width=0.2,
vertex.label=paste(names(V(stock_graph)), "-", V(stock_graph)$sector), # Asignar nombres que incluyan el sector
vertex.color = V(stock_graph)$color,
brush=TRUE, # resalte del nodo con un click
showLabels=TRUE, # mostrar nombres con hover
edge.alpha=0.6, # transparencia de las conexiones
bg="black", # background
main="Stocks graph using values with Graphical Lasso")
# save graph
graph_filename <- paste0("values_graph_rho_", rho, "_with_sectors.html")
saveWidget(stock_graph_js, file=graph_filename)
# open in browser
browseURL(graph_filename)